The role of spatial averaging in the precision of gene expression patterns 
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During embryonic development, differentiating cells respond via gene expression to positional cues 
from morphogen gradients. While gene expression is often highly erratic, embryonic development is 
precise. We show by theory and simulations that diffusion of the expressed protein can enhance the 
precision of its expression domain. While diffusion lessens the sharpness of the expression boundary, 
it also reduces super-Poissonian noise by washing out bursts of gene expression. Balancing these 
effects yields an optimal diffusion constant maximizing the precision of the expression domain. 



Embryonic development is driven by orderly, spatial 
patterns of gene expression that assign each cell in the 
embryo its particular fate. Experiments in recent years 
have vividly demonstrated that gene expression is of- 
ten highly stochastic Yet, embryonic development 
is exceedingly precise. A vivid example is the Bicoid- 
Hunchback system in the early Drosophila embryo, which 
has become a paradigm for understanding the formation 
of spatial gene-expression patterns. Shortly after fertil- 
ization, the morphogen protein Bicoid (Bed) forms an 
exponential concentration gradient along the anterior- 
posterior axis of the embryo, which provides positional 
information to the differentiating nuclei. One of the tar- 
get genes of Bed is hunchback (hb), which is expressed in 
the anterior half of the embryo. The posterior boundary 
of the hb expression domain is very sharp: by cell cycle 
13, the position of the boundary varies only by about 
one nuclear spacing [1, [H, This precision is higher 
than the best achievable precision for a time-averaging 
based read-out mechanism of the Bed gradient In- 
triguingly, the study of Gregor et al. revealed that the 
Hb concentrations in neighboring nuclei exhibit spatial 
correlations and the authors suggest that this implies a 
form of spatial averaging enhancing the precision of the 
posterior Hb boundary [3j. However, the mechanism for 
spatial averaging remained unclear. 

In this manuscript, we analytically and numerically 
study the Bcd-Hb system. Our analysis reveals a simple, 
yet powerful mechanism for spatial averaging, which is 
based on the diffusion of Hb itself. We show analytically 
that Hb diffusion between neighboring nuclei reduces the 
super-Poissonian part of the noise in its concentration, 
with a factor which depends on the diffusion length of 
Hb and the dimensionality of the system. In essence, dif- 
fusion reduces the noise by washing out bursts in gene 
expression. This mechanism is generic, and applies not 
only to any developmental system, but also to any bio- 
chemical network in general. For example, if a signaling 
protein is activated at one end of the cell and then has to 
diffuse to another place to activate another system, e.g. 
the messenger CheY in bacterial chemotaxis, then our 
results show that the non-Poissonian noise in the activa- 



tion of the signaling protein is washed out by diffusion; 
for this reason it may be beneficial to spatially separate 
the in- and output of a signaling pathway. Our analysis 
also reveals that while Hb diffusion reduces the noise.it 
also lessens the steepness of its expression boundary [3J . 
The interplay between these two antagonistic effects leads 
to an optimal diffusion constant of D ~ 0.1 /im 2 / s that 
maximizes the precision of the hb expression domain. 

We consider a minimal model of the Bcd-Hb system, 
which contains the necessary ingredients for understand- 
ing the mechanism of spatial averaging. This model in- 
cludes the stochastic and cooperative activation of Hb 
by Bed in each nucleus, and the diffusion of Hb between 
neighboring nuclei. To model the cooperative binding of 
Bed to the hb promoter, we assume that Bed proteins can 
bind sequentially to five binding sites on the hb promoter. 
When all five sites are occupied, the promoter is active 
(activity n = 1) and Hb proteins are produced stochas- 
tically with rate /3; otherwise the promoter is inactive 
(activity n = 0). Hb proteins are degraded stochasti- 
cally with rate fx. To obtain a lower bound on the pre- 
cision of the hb expression domain, we assume that Bed 
binds to the promoter with a diffusion limited on-rate 
k on = 4TraD B /V = 24<xD B /d 3 = 8.4 x 10 _5 s -1 where 
a ~ 3nm is the dimension of a binding site, d ~ 6.5 fmx 
the diameter of a nucleus and Db — 0.32 /im 2 / s the 
Bed diffusion constant Since the on-rate is assumed 
to be diffusion-limited, cooperativity of Hb activation is 
tuned via the Bed off-rate k Q s(j) — ajW which decreases 
with increasing number j of promoter-bound Bed. Using 
a = 12.3 s _1 and b = 6 the average promoter activity as 
function of the number of Bed molecules B in the nucleus 
approximately follows the Hill-function 



(n(B))~B h /(B h + K^), 



(1) 



with h = 5 and K eq = 690. The Hill-coefficient h 
was inferred from the relation between Hb and Bed 
concentration, and the activation threshold -Keq, where 
(n(.Keq)) = 0.5, is the average number of Bed molecules, 
(B), in a nucleus at the Hb domain boundary 0]. 

To describe the formation of the spatial Hb pattern, 
we place N = 64 x 64 = 4096 nuclei on a square lattice 
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with spacing t = 8.5 /zm [3[. Hb diffuses over the lat- 
tice with diffusion constant D. With reflecting bound- 
ary conditions in one and periodic in the other direction, 
our model is a cylinder which is symmetric around the 
anterior-posterior x-axis. Differential activation of Hb is 
induced by the Bed profile 



B(x) = K eq e 



-(x-x t )/>.t 



(2) 



which decreases exponentially with the distance x from 
the anterior pole at x = 0. The decay length is Ab = 
119.5 fj,m Q and the threshold position for Hb activation 
is set to x t = 0.5L; L — 64£ — 544 /im is the length 
of the embryo. Diffusion of Bed between nuclei induces 
fluctuations in the Bed copy number on the time-scale 
Td = ^ 2 /(4Db) — 10 2 s. Because Td is much smaller than 
the time-scale for promoter binding, ~ 10 s, Bed 
copy number fluctuations are efficiently averaged out by 
slow binding of Bed to the promoter. We therefore as- 
sume that the total number of Bed in a nucleus is con- 
stant and given by Eq. ([2]). The state of the system is 
described by the number of promoter bound Bed and 
the number of Hb in every nucleus. The state changes 
by reactions inside a nucleus and diffusion of Hb between 
nuclei. The dynamics is described by a reaction-diffusion 
master equation which we solve numerically [f| . 







<B>/Bo 

<n> 

---- <H>/Ho 


>file H/Ho 


1 

0.8 












O. 


0.6 




\ <m>:. 






ineou 


0.4 






(a) 






nstant; 


0.2 
n 


(b) 





anterior-to-posterior position x/L 



anterior-to-posterior position x/L 



FIG. 1: (a) Normalized average number of free Bed, {B)/Bo, 
average promoter activity, (n), and normalized average num- 
ber of Hb, (H)/Hq, in a nucleus as function of the position 
x/L on the anterior-posterior axis. The Hb decay rate is 
H = 1.2 X 10" 4 s _1 , the production rate p = 500^ = 0.06 s- 1 
and the Hb diffusion constant D — 0.32 fim 2 / s. The off-rate 
of the active promoter state is fc ff(5) ~ 1.6 x 10~ 3 s -1 . Nor- 
malization constants are Bo = 6480 and Ho = 500. (b) In- 
stantaneous normalized Hb profiles, H/Ho at different times 
in a steady-state trajectory for the same parameters as in (a). 
All simulations are for the full, two-dimensional system. 

Fig.[T^ shows simulation results for the average number 
of free Bed molecules, (B), the average promoter activ- 
ity, (n), and the average number of Hb molecules, (H) 7 
as function of the position x/L. Without Hb diffusion 
(H(x)) is proportional to (n(x)). For a finite Hb dif- 
fusion constant, however, the shape of the Hb profile is 
determined not only by (n(x)), but also by the diffusion 
length of Hb, A = y/ 1 D / ' \i — with increasing A the profile 
becomes less steep. Fig.QJi shows instantaneous Hb pro- 
files at different times. Fluctuations of promoter activity 
induce large fluctuations of H . These in turn lead to an 



uncertainty Ax in the position at which H crosses the 
threshold (H(x t )) = 0.5/3/m- 
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FIG. 2: (a) Simulation results for the Hb boundary width 
Aa; as a function of the Hb lifetime r = fi -1 for a constant 
production rate P = 0.06 s _1 and different diffusion constants 
D. The solid line scales as t~" 1//2 . (b) Simulation results for 
function of D and for P and /j, as in Fig.[T] The 
simulation results are compared to Eq. © using Eqs. (JS| and 
J5|. All simulations are for the full, two-dimensional system. 



We now consider the Hb boundary width when the Hb 
diffusion constant D is zero. Fig. [2^ shows that this de- 
creases with increasing lifetime of the Hb protein. To un- 
derstand this behavior, we note that the boundary width 
is, to a good approximation, given by 



Ax 



\{H{x t ))'\ 



(3) 



where a(xt) is the standard deviation of the Hb copy 
number and \{H(xt))'\ is the magnitude of the Hb gra- 
dient at the boundary position xt [1, 0| • When the Hb 
lifetime r is much longer than the time-scale fc ff(5) _1 of 
the promoter state fluctuations, then the noise in Hb copy 
number becomes Poissonian and the variance is given by 
the mean: a 2 (x t ) — (H{x t )) — (n(x t ))(3/n = 0.5/3r. On 
the other hand, the steepness of the boundary, |(ff(x t )}'|, 
increases linearly with r when D = and the synthe- 
sis rate ft is kept constant. Eq. © thus predicts that 
the boundary width Ax decreases as r -1 / 2 for large 
t 3> fc ff(5) _1 , which is indeed observed in the simu- 
lations (see Fig. [5^). Hence, the hb expression boundary 
could be made arbitrarily precise if the lifetime of the Hb 
protein could be increased indefinitely. 

In practice the averaging time cannot be made arbi- 
trarily long Q. Ultimately, it is limited by the nuclear 
division time r n d, which at cell cycle 13 is roughly 30 
minutes. Fig. [2^, shows that for r ~ r n d and vanishing 
Hb diffusion constant D = 0, the boundary width would 
be Ax ~ 4 nuclei, which is larger than the precision 
measured experimentally Q ■ The limited averaging time 
thus puts strong constraints on the precision that can be 
achieved via the mechanism of time-averaging alone. 

Fig. [2^ reveals, however, that the precision of the hb 
expression domain can be enhanced significantly by in- 
creasing the Hb diffusion constant to a finite value. This 
may seem surprising, since the steepness of the bound- 
ary, |(iJ(x t ))'|, decreases with increasing Hb diffusion 
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constant (Fig.QJ and this — as Eq. ([3]) shows — tends to 
increase the boundary width. However, increasing the 
Hb diffusion constant also reduces the noise in the Hb 
copy number. This is the mechanism of spatial averag- 
ing, which we now study analytically. 
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FIG. 3: (a) Hb is produced in nuclei and diffuses away from 
these Hb point sources. The stationary number of Hb in a 
given nucleus (shaded area) is the sum over contributions from 
all nuclei. The size of each contribution reduces exponentially 
with the distance from its source, (b) Slope (H(xt))' of the 
Hb profile at the boundary x t as function of the Hb diffusion 
length A = ^jD/yt. All lengths are in units of the lattice 
constant £ = 8.5 /jm. Analytical results from Eq. <(Sj (dashed 
line) are compared to simulation results (symbols) for the full, 
two-dimensional system. 

To elucidate the mechanism of spatial averaging, we 
first analyze how the steepness of the Hb boundary, 
\(H(xt))'\, depends on the Hb diffusion length A, and 
then how the noise in the Hb copy number at the bound- 
ary, <j(x t ), depends on A (see Eq. ([3])). To calculate the 
Hb profile, we note that a nucleus at a position Xi pro- 
duces Hb with an average rate f3(n(xi)) (Fig. [5^). In one 
dimension, the steady state number of Hb molecules in a 
nucleus of size I at x is approximately 



(H(x))t 



/3(n(xj))//i lx _ XiUX 
2{\/t) 



(4) 



The total Hb profile is the sum of (H(x)) ps over all point 
sources Xi = 0, £, 2£, . . . . To do the sum, we approximate 
the activity (n(xj)) from Eqs. |T]) and ([2]) by a piecewise 
linear function which behaves as {n(xi)} = 0.5 + m(xt — 
Xi) around the boundary for xt — 1 /2m < Xi < xt + I /2m 
and which is unity below and zero above the boundary. 
The slope — m = — /i/4Ab is the derivative of (n(xi)) at 
Xi — x t . At the boundary Xt the slope of {H(x)) is 



{H(xt))' = m^l-e- 
A* 



sinh 



1 



2mA 



(5) 



It increases linearly with the number of Hb molecules /3//J, 
and decreases as the Hb diffusion length increases. Due 
to the rotational symmetry, the slope in the full, two- 
dimensional system is also given by Eq. ([5]). Fig.[3j3 plots 
Eq. |5]) as function of A and compares the approximation 
with two-dimensional simulation results. It is seen that 
the agreement is very good. 
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FIG. 4: (a) Variance a^ B (x) of Hb as function of the dis- 
tance x from a point source at x = in one-dimensional 
space for increasing Hb diffusion length A. Eq. ((6)l is com- 
pared to simulation results for fi = 1.2 x 10 _2 s _1 and 
(3 — 500/i = 6s _1 , i.e. for fast Hb dynamics. While a pro- 
moter is active, /3/fc ff (5) ~ 3800 Hb molecules are produced 
on average. The Hb diffusion length is varied through the dif- 
fusion coefficient D. The inset shows the Hb noise strength 
<7p a {x '/ 'A) / '{H(x I 'A)) ps as function of distance normalized by 
the diffusion length, (b) Non-Poissonian part of the Hb noise 
strength, a 2 / (H) — 1, for a uniform system with N = 32 x 32 
nuclei on a square lattice as function of X/i. Decay and pro- 
duction rate of Hb are constant, while D is varied; for other 
parameter values, see Fig.[T] The average promoter activity 
is (n) = 0.5. Simulation results for the uniform system (cir- 
cles) are compared to Eq. (JSJ with the two-dimensional value 
AT (A) = 8(A/£) 2 (solid line). The value a$/(H) - 1 = 82.4 (see 
Eq. in an isolated volume has been calculated numerically. 
These results are compared to the Hb noise strength for the 
full, two-dimensional system (see Fig.[T]) in a nucleus at the 
boundary with (n) = 0.5 (triangles). 



Next, we study the effect of Hb diffusion on the vari- 
ance in the Hb copy number, a 2 (see Eq. ([3|)). For clarity, 
we first do the calculations for one-dimensional diffusion 
and extend our approximation to general dimensions af- 
terwards. To compute tr 2 , we exploit the observation 
made above that each nucleus acts as a point source of 
Hb (Fig. [3^), and that the expression of Hb in each nu- 
cleus is an independent stochastic process. We denote 
the variance of the Hb copy number in a nucleus of size 
t at a distance x from the point source as a^ s {x). In a 
uniform space, the total Hb variance in a nucleus is the 
sum a 2 — a 2 s (xi) over point sources at all Xi. To cal- 
culate <jp S (a;), we assume that the Hb dynamics is fast on 
the time-scale of promoter switching, such that the Hb 
concentration switches between zero and (H(x)) ps /(n), 
with (H(x)) ps given by Eq. ((4J. In this limit, 



a 2 ps (x) « (H(x)) ps + (H(x))l 



(n) : 



(6) 



where the second term is due to promoter switching with 
variance (1 — n)(n), and the first term describes the Pois- 
son noise coming from the production, diffusion and de- 
cay of Hb when the promoter is active. 

Fig. UK shows Up S {x) for different values of the diffu- 
sion length A. For small A, <J 2 s (x) is large for small x, 
but reduces quickly with increasing x. Hence, the vari- 
ance at a given nucleus is determined by few nuclei in 
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the immediate neighborhood but each individual contri- 
bution is large. For increasing A the number of nuclei 
contributing to the variance increases, but the individual 
contributions of nearby nuclei are smaller; this is because 
diffusion washes out the bursts of Hb production at the 
source. To see which of the two opposing effects domi- 
nates, we sum Eq. (|6]) over X4 to obtain the total variance: 



(H) 



(H) 2 (l-n) 
4(X/£) (n) 



(7) 



Compared to the variance in an isolated nucleus, the 
overall effect of diffusion in a one-dimensional, uniform 
space is a reduction of the non-Poissonian noise by 
(4A/£) _1 . The Poissonian part cannot be reduced be- 
cause diffusion itself is Poissonian. Indeed, spatial aver- 
aging can reduce the effect of noise in gene expression, 
but only if this is super-Poissonian. Simulations also al- 
low us to compute a 2 s (x) for slower Hb dynamics. In this 
case, the non-Poissonian part of the variance is smaller. 
However, the qualitative dependence of x remains the 
same: the non-Poissonian part is washed out by diffu- 
sion over a distance proportional to A (inset Fig. 2^,). 

The finite range of non-Poissonian fluctuations im- 
plies the following simple model. N nuclei contribute 
equally to a well-stirred reservoir containing M Hb 
molecules, with average (M) = N(H) and variance 
a 2 M = N(Jq, where Uq is the noise of H in a nucleus 
if there was no Hb diffusion. Each nucleus then sam- 
ples Hb molecules in a binomial fashion from this reser- 
voir. The noise in each nucleus is formally given by 
a 2 =E[V[H\M}]+ V[E[H\M}}, where E[X] denotes the 
expectation and V[X] the variance of X. The first term 
describes the noise in the sampling process, and is given 
by E[V[H\M}\ = (H)(l - N^ 1 ). The second term de- 
scribes the variance in H due to fluctuations in the reser- 
voir, and is given by F[i?[iJ|M]] = cr 2 /N. Hence, 



a 2 = (H) 



°l ~ (H) 
N 



(8) 



This simple expression elucidates that it is the super- 
Poissonian part of the noise, a 2 — (H), which is reduced 
by spatial averaging. Since N scales with the dimension- 
ality d as N ~ X d , this reduction becomes more efficient 
in higher dimensions. Interestingly, this expression also 
reveals that if hb expression would be sM&-Poissonian Q , 
diffusion would increase the noise. 

Fig. [4)3 shows the non-Poissonian part of the variance 
in a uniform, two-dimensional system as function of the 
diffusion length; the exact, numerical result is compared 
to the prediction of Eq. ([8]) with N chosen to be the num- 
ber of nuclei within a distance 2A of a given nucleus along 
the edges of a square lattice, N(X) = 8(X/£) 2 . For suffi- 
ciently large A for which N(X) > 1, the approximation is 
excellent. Fig.fJJ) also compares the results for the uni- 
form system to those of the full, non-uniform system in 



which a Bed gradient activates Hb. Clearly, the well- 
stirred approximation describes the variance at the Hb 
boundary very well; only for A > 51 is the variance at 
the boundary significantly influenced by nuclei far away 
from the boundary with smaller variance. 

Finally, we can combine Eqs. ([5|) and (J8]) with Eq. (J3j> 
to predict the boundary width Ax. Fig. (2b shows that 
the prediction agrees very well with the simulations. This 
figure also shows that the two antagonistic effects of Hb 
diffusion — reducing the slope but also the variance of the 
Hb concentration — lead to a diffusion constant that opti- 
mizes the boundary width for a fixed Hb lifetime. Inter- 
estingly, this minimal width is less than one nuclear spac- 
ing, as found experimentally Q. This suggests that for 
achieving the necessary precision, mechanisms based on 
multiple gradients or interactions between multiple Bed 
targets are not required, although these may provide ro- 
bustness against embryo-to-embryo variations 0, [ToL 11 1 . 
Since the exact mechanism for Hb transport is unknown, 
it is tempting to speculate that the precision achieved via 
spatial averaging could be increased further by separating 
the two effects of Hb diffusion using an anisotropic trans- 
port of Hb: slow transport along the anterior-posterior 
axis would allow for steep spatial profiles, while fast dif- 
fusion along the perpendicular direction would allow for 
effective spatial averaging. The benefit of anisotropic 
transport could imply the existence of active mechanisms 
for Hb transport. 

Spatial averaging can only be beneficial when the noise 
in Hb production has a super-Poissonian component. 
The observation of spatial correlations in Hb by Gregor et 
al. indicates that this is the case Q. The question arises 
whether bursts in gene expression are inevitable and spa- 
tial averaging a prerequisite for achieving a precise hb 
expression domain. The Hb lifetime is limited by the 
nuclear division time, while the time-scale for promoter 
state fluctuations is limited by the diffusion of Bed. The 
magnitude of the bursts could be reduced by reducing 
the promoter strength, but this would lower the Hb copy 
number and hence the steepness of the Hb boundary; the 
net result would, in fact, be a decrease of its precision. 
In a model with explicit hb mRNA, the effect of pro- 
moter state fluctuations could be alleviated by reducing 
the transcription rate. However, to achieve a sufficiently 
steep and precise Hb boundary with a low number of hb 
mRNA, the translation rate had to be increased such that 
Hb production becomes super-Poissonian again. It thus 
appears that spatial averaging is a fundamental mecha- 
nism for generating precise patterns of gene expression. 
Clearly, it will be of interest to study how other modes 
of transport, e.g. sub-diffusion or active transport, affect 
the mechanism of spatial averaging. 
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